Dissipation and Ultrastrong Coupling in Circuit QED 
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Cavity and circuit QED study light-matter interaction at its most fundamental level. Yet, this 
interaction is most often neglected when considering the coupling of this system with an environment. 
In this paper, we show how this simplification, which leads to the standard quantum optics master 
equation, is at the root of unphysical effects. Including qubit relaxation and dephasing, and cavity 
relaxation, we derive a master equation that takes into account the qubit-resonator coupling. Special 
attention is given to the ultrastrong coupling regime, where the failure of the quantum optical 
master equation is manifest. In this situation, our model predicts an asymmetry in the vacuum 
Rabi splitting that could be used to probe dephasing noise at unexplored frequencies. We also show 
how fluctuations in the qubit frequency can cause sideband transitions, squeezing, and Casimir-like 
photon generation. 

PACS numbers: 42.50.Pq, 03.65.Yz, 42.50.Lc 



I. INTRODUCTION 

Elementary quantum mechanics teaches that a closed 
physical system always evolves in a reversible manner. 
However, control and readout imply a coupling of the 
quantum system to the outside world, making it sub- 
ject to relaxation and decoherence. These irreversible 
dynamics, well understood from a theoretical viewpoint, 
have also been experimentally tested. In cavity QED lj, 
macroscopic superpositions of quantum states of light 
have been built and their destruction due to their in- 
teraction with a reservoir has been observed [2]. Using 
repeated QND measurements, the birth and death of sin- 
gle photons in a cavity has been studied [3J dj . Circuit 
QED, a solid-state realization of cavity QED, also offers a 
detailed understanding of relaxation and dephasing phe- 
nomena [3 |S] . Spontaneous emission of a qubit in a res- 
onator has been characterized with respect to the influ- 
ence of far off- resonant modes [7]. Moreover, the impact 
of measurement on qubit dephasing processes is well un- 
derstood [HI W\ , for instance in the cases of dispersive [TO] 
and bifurcation [llj read-out. 

Though both circuit and cavity QED allow to study 
dissipation, solid-state devices allow much stronger light- 
matter interaction rates. For example, current-current 
coupling of a flux qubit to a Josephson junction in a res- 
onator can boost the strength up to the order of the res- 
onator and qubit frequencies |12j . breaking the rotating- 
wave approximation (RWA). This ultrastrong coupling 
regime has been achieved experimentally with coplanar 
waveguides [T3] and lumped LC resonators Q3] . In paral- 
lel to these experimental efforts, dynamics of pure states 
have been theoretically studied [15H18] , A rigorous model 
based on the Bloch-Redfield formalism which describes 
photon losses has also been proposed by Hausinger and 
Grifoni [TO] . Finally, a non-Markovian model of dissipa- 
tion has been used to predict the emission spectrum of an 
atom-resonator system in which the ultrastrong coupling 
strength is modulated over time [201 HI] an d to study 



the sensitivity of the system to noise in the qubit fre- 
quency [22] . In this paper, we give a complete descrip- 
tion of dissipation including qubit relaxation and pure 
dephasing in the ultrastrong coupling regime, focusing 
on the standard case where the baths can be treated as 
Markovian. 

Qubit-resonator coupling is at the heart of the prob- 
lem with dissipation in the ultrastrong coupling regime. 
When the coupling between these two subsystems is 
small, interactions with the environment are treated sep- 
arately for the qubit and the oscillator pQ. However, 
when the atom-field interaction increases up to the break- 
down of the RWA, this approach leads to unphysical pre- 
dictions. For example, and as will be illustrated later 
in Fig. [T] relaxation baths bring the system out of its 
ground state even at T = 0. Furthermore, in the pres- 
ence of a strong qubit-resonator coupling, transitions at 
widely separated frequencies appear, breaking down the 
standard white noise approximation. To avoid such an- 
noyances, the qubit-resonator coupling and colored baths 
must be included in the treatment of dissipation. 

The outline of this work is as follows. First, we present 
the system Hamiltonian, along with a perturbative ap- 
proach that can diagonalize it approximately. Next, we 
discuss the treatment of dissipation. In Section III A 



we explain the standard approach to describe dissipation 
in the Jaynes-Cummings regime. We then show issues 
arising from the use of this technique in the ultrastrong 



regime, and devote Section III B to the presentation of a 
Lindbladian master equation that solves them. In Sec- 
tion |IV| we describe the implication of these results in the 
strong coupling regime. Finally, in Section [Vl we study 
physical consequences of the model obtained here in the 
Bloch-Siegert regime, for which counter-rotating terms 
can be treated in a perturbative fashion. We first show 
how the vacuum Rabi splitting spectrum is affected by 
non-RWA terms in the Hamiltonian and by the shape of 
the noise spectrum. We also introduce a potential tech- 
nique to exploit these effects to study noise. We then 



present how the master equation helps to understand an 
analog of the time-dependent Casimir effect coming from 
pure qubit dephasing [25] , 



II. HAMILTONIANS DESCRIBING THE 
QUBIT-RESONATOR SYSTEM 

The Rabi Hamiltonian, describing the interaction of a 
two-level atom with a single electromagnetic mode of a 
resonator, takes the form [h = 1) [23] 



Hr = u> r a/ a 



-er 2 + gXa x , 



(1) 



where uj a is the qubit splitting, ui r the resonator fre- 
quency, g the coupling strength, and X = a) +a. In most 
experimental situations, g <C ui a ,uj r and the rotating- 
wave approximation (RWA) can safely be made. This 
amounts to dropping the fast-oscillating, or counter- 
rotating, terms Iqk = a<j-+aSo+ from Hr. This approx- 
imation leads to the Jaynes-Cummings Hamiltonian |25j 

Hjc — Lo r a)a + -^&z + g(a<?+ + aV_). (2) 

In opposition to the Rabi Hamiltonian, here the total 
number of quanta N q = (1 + cr z )/2 + a^a is a good 
quantum number, allowing exact diagonalization of Hjc- 
The system enters the ultrastrong coupling regime when 
g is so large with respect to uj a , uj r that Iqr leads to 
experimentally observable consequences and the RWA 
cannot be safely made |T3l [14] . In this situation, since 
[N q , Icr] ^ 0, the total number of excitations is not pre- 
served, even though its parity is [15 j . As a result, even in 
the ground state, the expected mean number of resonator 
and qubit excitations is non-zero. 

Although the analytical spectrum of Hr has recently 
been found by Braak [26], it is defined in terms of the 
power series of a transcendental function. An approxi- 
mate, but more simple form, can be found in the interme- 
diate regime where g is small with respect to £ = ui a -\-u r , 
with the system still being in the ultrastrong coupling 
regime. This will be referred to as the Bloch-Siegert 
regime. This is done using the unitary transforma- 
tion omiiT] 

U = exp {A(acr_ - aV+) + £(a 2 - a t2 )a z } , (3) 

where A = gfS, and £ — gA/2ui r together with the 
Campbell-Baker-Hausdorff relation 



c- x Hc x = H 



[H,X] + ±[[H,X},X} + ... (4) 



To second order in A, this yields the Bloch-Siegert Hamil- 
tonian 

U^H R U ~ H BS = (« P + ficr z )a^a + ^-a z + gl+, (5) 



where /+ = aer + -t-aV_, u) q = u> a +n, and /i = ,g 2 /S. This 
Hamiltonian is similar to the Jaynes-Cummings Hamil- 
tonian, but contains Bloch-Siegert shifts fj, on qubit and 
resonator frequencies. 

Since Hqs is block-diagonal, its eigenstates can be 
found exactly to be 

\n, +) = — sin6*„|e,n — 1) + cos 9 n \g,n) (6) 

\n,—)= cos9 n \e,n- I) +sm9 n \g,n), (7) 

with the Bloch-Siegert mixing angle 



arctan 



Ag s - v /(AgS)2+4^ n 
2gy/n 



(8) 



and where A^ s — uj a — uj r + 2pm. To second order 
in A = .g/S, the excited eigenstates \n, ±) of the Rabi 
Hamiltonian in the bare basis are then given by 



\n,±) = U\n,±), 

while the ground state takes the form 

A 2 
2 



|. 9 0> = U\g0) 



1 



(9) 



\ 9 0)-A\el)+^V2\g2). (10) 



As mentioned before, the ground state is no longer the 
simple Hjc ground state \g0), but now contains qubit- 
resonator excitations. 

Unitary transformation Eq. ([3]) deserves further at- 
tention. With the replacement a± — > cr*', the term 
proportional to A generates a displacement of the res- 
onator field. Moreover, the term proportional to £ gener- 
ates squeezing of the field, with a qubit-state dependent 
squeezing parameter £. We can thus expect the qubit- 
resonator state to display the properties of displaced- 
squeezed states, both transformations being qubit-state 
dependent. With £ = gA/2u r , squeezing is expected to 
be larger for uj r <C Lo a [28 . 



III. MASTER EQUATIONS 

In this section, we introduce dissipation following two 
approaches. First, we follow the standard approach 
where the qubit-resonator coupling g is ignored when ob- 
taining the master equation [T| . This results in the stan- 
dard quantum optics master equation |25) . We then con- 
sider an approach taking into account the non-negligible 
qubit-resonator coupling. In both cases, the qubit and 
the resonator are assumed to be weakly coupled to a bath 
of harmonic oscillators, with free Hamiltonian 



H B = J2"ibl b i 



(11) 



where 6; , b\ are ladder operators for bath mode I with 
frequency vi and system-bath coupling 

H SB =Y,Mc + ^)(bi+b]), (12) 



with ati a coupling strength to bath mode I. For the 
qubit c —¥ er_, while for the resonator c — > a. In the 
standard approach, this will correspond to qubit and res- 
onator damping, respectively. Finally, dephasing is mod- 
eled classically as 



Hdep = f(t)cr z , 



(13) 



where f(t) is a random function of t with zero mean value. 
A quantum model for dephasing leads to similar results 
and is presented in Appendix |B 2| for completeness. As 
will be seen in Section |VJ though the master equations 
obtained in the quantum and the classical cases have the 
same form, asymmetric noise spectral densities are al- 
lowed in the quantum model, yielding different predic- 
tions in the ultrastrong coupling regime. 



A. Standard master equation 

The standard approach is to assume that qubit and 
resonator are independent when obtaining the dissipa- 
tive part of the master equation. The coupling is then 
reintroduced in an ad-hoc fashion in the Hamiltonian. 
This leads to the standard master equation 



dp 
df 



-i [H, p] + £ s tdP, 



where, at T = 0, 



£ st d- = n,V[a] ■ +7i£>[ct_] 



l<t> 



V[a z }- 



(14) 



(15) 



with V[0] P = \ (20 pO* - pO^O - O^Op). Here, k is 
the photon leakage rate for the resonator, 71 the qubit 
relaxation rate and 7^ the qubit pure dephasing rate. 
This expression is obtained in the Markov approximation 
which assumes the spectral density of all three baths to be 
white. In other words, the environment-system coupling 
is evaluated at the relevant frequency {u r for k, ui a for 
71 and uj — > for 7^) and then assumed to have support 
at all frequencies. 

For g/T, small enough for the RWA to be safely per- 
formed, this expression while not rigorous 29, [3U] can be 
used to accurately describe many cavity QED and cir- 
cuit QED experiments [TJ [BJ. Indeed, the terms T>[a\- 
and T>[aJ\- in Eq. ( [15] ) correctly tend to bring the sys- 
tem to the ground state |g0) of the Jaynes-Cummings 
Hamiltonian. 

In the ultra-strong coupling regime however, \gQ) is 
no longer the ground state and Eq. ( 14 ) will bring the 



ultrastrongly coupled qubit-resonator system outside of 
its true ground state |<?0). Therefore, even at T — 0, in 
which case no energy should be added to the system, re- 
laxation will generate photons in excess to those already 
present in the ground state. These additional excitations 
are plotted in Fig. [I] as a function of g as the black line. 
This curve closely follows the behavior of the error one 
makes by approximating the Rabi ground state by the 




1000 1500 
9 (MHz) 

FIG. 1. (Color online) Excess in the mean photon number 
due to relaxation in the steady state of the ultrastrong qubit- 
resonator system. Initially, the system is in its true gro und 
state \g0), but, under the standard master equation (14 1, re- 
laxation unphysically excites the system even at T = 0. The 
black line, which corresponds to the left axis, represents the 
number of additional photons introduced in steady state by 
dissipation. The red dots, associated to the right axis, desig- 
nate one minus the fidelity of the Rabi ground state \g0) to the 
vacuum state |<?0). The parameters are ui a /2ii — uj r /2ir — 6 
GHz, n/2ir = 7i/2-7r = 0.1 MHz and no pure dephasing. In- 
set: mean photon number as a function of time for the system 
starting in its ground state with g/2n = 2 GHz. In both the 
main plot and the inset, the blue dashed line indicates results 
for the fidelity and the photon number as obtained with the 
master equation presented in Sec. |IIIB| 



vacuum state, which is represented by the red dots. It is 
important to emphasize that these results are obtained 
for an undriven system evolving simply under dissipative 
dynamics. We also note that, in practice, preparing \g0) 
can be extremely challenging, requiring for example tun- 
ing of the coupling constant g in a time scale <C 1/5, this 
being typically in the sub-nanosecond range. This is why 
the system is initialized in \g0) in Fig. Ill 



Moreover, because it assumes white noise, Eq. ( 15 ) in 



correctly describes Purcell decay, which can be the factor 
limiting coherence in superconducting qubits [7 . Addi- 
tionally, Purcell decay is probing the resonator bath at 
the qubit transition frequency, something which is miss- 
ing from the above description but is central to the exper- 
iments reported in Refs. [3TJ [35] For the same reason and 
as discussed below, it also incorrectly describes dressed- 
dephasing [331I34] . 



B. Master equation in the dressed picture 

We now take into account qubit-resonator coupling 
when deriving the master equation. In this case, we can- 
not assign a unique dissipation channel to each bath men- 
tioned above. Indeed, rather than transitions between 
eigenstates of the free Hamiltonian Ho = oj r a'a+uj a a z /2, 
coupling to the baths leads to transitions between the 
qubit-resonator entangled eigenstates {\n, ±), \g0)}. To 
simplify the notation, these states will be denoted below 
as \j), j increasing with energy. These states can be ap- 



proximated analytically as explained in Sec. [IT] or found 
exactly numerically ^7\ I35H38] . 

To obtain a master equation that takes into account 
the coupling <?, we first move to the frame that diagonal- 
izes the Rabi Hamiltonian for both the system and the 
system-bath Hamiltonians. Neglecting high-frequency 
terms, the resulting expressions involve transitions \j) o 
|fc) between eigenstates at a rate which depends on the 
noise spectral density at frequency Akj = u>k — Wj. If 
their linewidth is small enough, these transitions can be 
treated as due to independent baths. As a result, these 
independant baths can each be treated in the Markov 
approximation [39 . As shown in Appendices [A] and [Bl 
this leads at T = 0, to the Lindbladian 



J l II II II " 



X and a r baths <x- bath 



|2+> 
12-) 



FIG. 2. (Color online) Transitions driven by noise. X and a x 
baths can only generate transitions between states of different 
parity. The a z bath can generate transitions between any pair 
of levels of same parity. 



£dr 



V 



E^IjXjI 



■ £ if ©[|j)<fc|]. 

j-k^j 



■ Yi (n k +ri r k )v[\j)(k\]. 

j,k>j 



(16) 



where \j) and \k) are eigenstates of the qubit-resonator 
system. Temperature dependence is taken into account 
in the Appendices but dropped here to simplify the dis- 
cussion. The first two terms in Eq. ( 16 ) are the contri- 



butions from the bath described by Eq. ( 13 ) that caused 
only dephasing in the standard master equation. Here, 
this a z bath causes dephasing in the eigenstate basis with 



*,= 



7^(0) 



(17) 



where 7^(w) is the rate corresponding to the dephasing 
noise spectral density at frequency uj and where 



Jk _ 



UWz\k). 



(18) 



Since a z is not diagonal in the eigenbasis, it also causes 
unwanted transitions at a rate 



,jk _ l4>( A kj) I ,7c I 



(19) 



This contribution will only be significant if the dephasing 
bath has spectral weight at the potentially high frequency 
Akj or if the qubit is operated away from a sweet-spot, 
in which case even low spectral weight can have a large 
impact [ID]. Finally, the last two terms of Eq. ( p~6[ ) are 
contributions from the resonator and qubit baths that 
caused energy relaxation in the quantum optical mas- 
ter equation. They now cause transitions between eigen- 
states at rates 



where 



ri k = K (A kj )x \x* k \ 

F 7 fc = 7 (A fej -)x|4f 



x° k = 01*1*) 

4 k = (j\a x \k). 



(20) 

(21) 

(22) 
(23) 



Here, k(u>) and 7(0;) are rates that are proportional to 
noise spectra, respectively for resonator and qubit envi- 
ronments. 

The dressed Lindbladian £dr ■ solves the problem stated 
in Section |IH A| Indeed, at T = 0, rather than excit- 
ing the system, dissipators accounting for relaxation in 
Eq. ( 16 1 lead to decay to the true ground state. This 
is illustrated by the dashed blue line in Fig. [I] More- 
over, it is interesting to point out that, in addition to 
the zero-frequency term responsible for pure dephasing, 
the noise along a z can stimulate transitions between the 
eigenstates \j), leading to dephasing-induced generation 
of photons and qubit excitations |4T] . This is related to 
the time-dependent Casimir effect, as will be discussed 
further in Section [V Bl 

Finally, Fig. [2] illustrates the allowed transitions given 
the symmetry of the Rabi Hamiltonian, and in particular 
given that it preserves the parity of the total number of 
excitation. As further explained in Appendix [A] for odd 
(parity-changing) transition matrices, such as relaxation- 
related operators X and a x , no decay is possible be- 
tween states of same parity. On the other hand, the even 
(parity-preserving) a z matrix associated with dephasing 
can generate transitions only between pairs of states of 
same parity. 



IV. STRONG COUPLING REGIME 



Before going to the ultra-strong coupling regime, in 
this section we focus on the simpler strong-coupling 
regime described by the Jaynes-Cummings Hamilto- 
nian. We first consider the dispersive regime, a situ- 
ation which is particularly useful for qubit readout in 
circuit QED, and then move to the Jaynes-Cummings 
Hamiltonian. The analysis done here is in the spirit of 
the dressed-dephasing model [53], but also encompasses 
qubit-resonator resonance [42] . 



A. The dispersive regime 

As already mentioned above, due to the white noise 
approximation, the description provided by the stan- 
dard master equation of Sec. |IH A| can break down 
even when the dispersive approximation is valid. In 
this regime, achieved when |A| = \uj a — u r \ ^> g, the 
Jaynes-Cummings Hamiltonian reduces to the dispersive 
Hamilonian 



-ffdisp = (w r +x°"z + C)a t a+- 



X 



Z + Qa z (a^a)\ (24) 



to fourth order in g. This diagonal Hamiltonian includes 
an effective qubit-resonator dispersive coupling x = 9 2 1 ^ 
and a small nonlinearity C, = <? 4 /A 3 which is usually ne- 
glected. The coupling implies that the dispersive eigen- 
states display some degree of mixing between qubit and 
resonato r. In deed , to s econd order in g, these eigenstates, 
denoted \gn) and \en), axe 



\e,n- 1) 



1 



9 2 n 

2A 2 



e,n - 1) 



gyn , 

A 



\9,n) 



\9,n) 



g 2 n 



9V n 



e, n — 1) 



(25) 
(26) 



Consequences of this mixing of qubit and resonator 
states are Purcell decay [43] and the dressed-dephasing 
model [33] EH]- Purcell decay is the relaxation of the 
qubit by photon emission out of the cavity. This is cap- 
tured here by the rate T^ evaluated between the dressed 
states |e,0) and \gffj yielding Tf >s° = (g/A) 2 n(uj a + x ). 
In this expression, the standard Lindbladian Eq. ( 15 ) 



rather evaluates cavity damping at the cavity frequency 
lj t . This difference is important in several circuit QED 
experiments [TI I3T1 132"] . 

The standard approach also does not capture cor- 
rectly dressed-dephasing discussed in Refs. [33j[34]. Es- 
sentially, dressed-dephasing captures how dephasing can 
produce relaxation because of the finite qubit-photon 
mixing in the dressed states Eq. ( 25 1 . Deriving the 



dressed-dephasing rate from Eq. (15 1 yields a result pro- 



portional to the spectrum of dephasing noise at zero fre- 
quency. Here and in Refs. [33JEH], we rather obtain the 
rate Tl involving the spectrum of dephasing noise at the 
qubit-resonator detuning frequency. Assuming a dephas- 
ing noise scaling as 1//, the difference between the two 
predictions can be quite large. In practice, this means 
that one must be careful in interpreting results of nu- 
merical simulations of the standard Lindbladian Eq. ( 14 ) 



as it can include unrealistically large amounts of qubit 
flipping induced by dephasing noise. It is worth point- 
ing out however that dressed-dephasing can be relevant 
experimentally in some circumstances |40) . 

In Refs. [33J EH], the nonlinear term proportional to 
C is neglected in the derivation of the master equation 
containing the dressed-dephasing contribution. This ap- 
proximation breaks down when £ > k, in which case the 
approach developed here is appropriate. However, when 



£ ~ k, the different resonator transitions are not well 
separated and the environment cannot be treated as in- 
dependent baths. The approach of Refs. [33J EH] should 
then be used. The validity of the results obtained here is 
further discussed in Appendix [C] 



B. Dissipation in the Jaynes-Cummings model 

In this section, we consider the situation where the 
dispersive approximation does not hold but the RWA is 
still valid. This is essentially generalizing the results of 
the dressed-dephasing model [33JEH]. Under the RWA, 
the ground state is simply \g0). Excited eigenstates |n±) 
are given by Eqs. ([6]) and (17]), with the mixing angle 0„ 
defined by Eq. Q and /i = such that A^ s = uj a — to r — 
A. 

We first consider the matrices X and a x , whose el- 
ements are involved in relaxation rates described by 
Eqs. (20) and (21). To keep the discussion simple, we 
limit ourselves to the subspace {\g0), |1— ), |1+)}- Com- 
plete results can be found in Appendix ID] Since the 



Jaynes-Cummings eigenstates have a well-defined exci- 
tation number, only transitions involving the loss or gain 
of one quantum are allowed, thus forbidding transitions 
between |1+) and |1— ). This yields 



Y = | Miitf, 1 (27) 

(28) 



With the eigenstates changing character with 9\ between 
mostly qubit-like or photon-like, the contribution of the 
two decay channels X and a x follows. This can be visual- 
ized geometrically as illustrated in Fig. [3k. In particular, 
when the qubit and the resonator are on resonance, their 
corresponding relaxation noises have exactly the same 
weight. For example, the matrix elements of X reduce 
to 



sin 9i cos 0i \ 


sin 0i 





cos 0i 


° J 


cos 0i 


— sin 0i 


cos 0i 





-sin 0i 






X 90: 



1,± 



±- 



1 



y , n ) +;n+l J + y*n,-;n+l,- 

Y*n,+;7i+l,— •yn,—;n+X,-\ 



(29) 

(Vn + Vn + 1) (30) 
(Vn - Vn + 1) , (31) 



which exactly leads to the master equation presented in 
Ref. [12] in the presence of resonator losses only. 

Under the RWA, a J z k can only be non-zero for states 
that involve the same total excitation number, i.e. are 
in the same Jaynes-Cummings doublet. The resulting 



a) JeO) 




FIG. 3. (Color online) Matrix elements under the RWA in the 
subspace {|<?0), |1— ), |1+)}- a) Relaxation matrix elements. 
Eigenstates are a mixture of qubit and resonator states, with 
an angle B\. The fraction of the relaxation rate that comes 
from the qubit or the resonator bath is determined by the pro- 
jection of the eigenstate on the qubit |e0) or resonator \gl) 
axis, b) Dephasing matrix elements. The dephasing Hamil- 
tonian rotates state vectors around the |e0) axis. Resulting 
vectors have a projection on the orthogonal eigenstate in the 
same doublet. This generates transitions between |1,+) and 
|l,-)if0i>O. 
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FIG. 4. (Color online) Vacuum Rabi splitting, a) Transition 
rates involved in the perturbative calculation in the three-level 
approximation, b) Schematic plot of Im(a) s as a function of 
aid- In general, the result is not symmetric. 



matrix elements are 



a 



7 so ;9 o 



7i±:n± 



a 



n=p;n± 



Tcos(26> n ) 
— 2cos6L sin( 



(32) 
(33) 
(34) 



Generalizing the dressed-dephasing model, the above for- 
mulae show that the dephasing bath induces transitions 
between states in the same JC doublet. As illustrated in 
Fig.[3j}, this happens only if there is some mixing between 
the qubit and the resonator. In particular, in resonance, 
a n±-,n± _ q anc j a n+-,n± _ ^ Then, dephasing processes 

for states that do not involve \g0) are entirely due to 
transitions within the doublets which are caused by de- 
phasing noise at the doublet splitting frequencies 2gy / n. 
Since these are very high frequencies and dephasing is 
often caused by a 1/f bath [44], transition rates within 
doublets are expected to be small. Therefore, in reso- 
nance, states that do not contain \g0) should be largely 
immune to pure phase-destroying processes. 



V. ULTRASTRONG COUPLING REGIME 

Here, we take the ratio g/Y, to be sufficient to break 
the RWA, but still much smaller than unity. In this situ- 
ation, the transition matrix elements given in Eqs. (18), 



(22), and (23) can be evaluated using the Bloch-Siegert 
eigenstates Eqs. ^ and Q, as done in Appendix |d| In 
this section, we use these results to study two distinctive 
phenomena occuring in the ultrastrong coupling regime: 

1) asymmetry of the vacuum Rabi splitting spectrum; 

2) sideband transitions and photon generation caused by 
qubit frequency modulations. 



A. Asymmetry of the vacuum Rabi splitting 
spectrum 

Vacuum Rabi splitting is observed by measuring trans- 
mission (Im (a) and/or Me (a)) of the resonator under 
weak rf excitation [B] [[5] . In the presence of the rf drive, 
the Hamiltonian becomes 



Hdrvn(t) = H R + eac l 



e a'e 



f -iu d t 



(35) 



with e the amplitude of the drive and uJd its frequency. 
Assuming j< S, we find in Appendix |F| under the three- 
level approximation that Im (a) in steady-state and at the 
Bloch-Siegert-shifted qubit-resonator resonance (A BS = 
0) is given by 



eri/2 



er 2 /2 



rf + (A B S + gf 



r2 + (A BS -.g) 2 ' (36) 



where 



r i = g (J- +70 + 70 



r 2 = 2 U 



(7+ +7^ +7^) 



(37) 
(38) 



The various rates entering these expressions are illus- 
trated in Fig. Ek) and can be found in Appendix IfJ As 
in the standard case, the transmission is composed of 
two Lorentzians separated by 2g Tj. However, here two 
distinct rates Ti and T 2 dictate the width and height 
of these peaks. As a result, the vacuum Rabi splitting 
spectrum can be asymmetric, even when the qubit and 
the resonator are in resonance. Asymmetry in the pres- 
ence of counter-rotating terms has also been pointed out 



considering only cavity decay in the Fourier transfom of 
(cr z (i)) [19] and qubit fluorescence [46] . 
Here, three situations can lead to asymmetry: 

1. Relaxation noise spectra are not equal at the 
frequencies corresponding to the two transitions 
|1±) — > \g0). This situation will be referred to as 
the non-white relaxation noise case. 

2. The pure qubit dephasing noise spectrum is not 
equal at frequencies Ai± i T . Since classical noise 
spectra are always symmetric in frequency |47j . 
we call this situation the quantum dephasing noise 
case. 

3. Keeping counter-rotating terms such that A 5/ 0, 
the matrix^elements_of X and o~ x for transitions 
|1±) — > |g0) and |1±) —t | 1=f) are not equal, as 
shown in appendix [D] This is the ultrastrong case. 

To characterize the asymmetry as the coupling g is 
increased, the noise spectra must be known. We now 
make some hypotheses on that noise and consider their 
consequences. We first isolate the influence of counter- 
rotating terms by choosing white relaxation noise and no 
pure dephasing. In these conditions, the second order 
terms in g cancel out and the asymmetry r/ = T± — T2 
increases linearly with g: 



A 



>h 



(k + 7i). 



(39) 



For the parameters realized in Ref. 13J, g/2ir — 
636 MHz, u r /2n = 5.357 GHz, and k/2-k = 3.7 MHz, 
and taking 7i/2tt = 0.1 MHz yields t)/2it ~ 0.11 MHz 
and in turn an asymmetry of ~ 6% in the transmission 
peak amplitudes. As a result, in the ultrastrong regime, 
the height of the transmission peaks in a vacuum Rabi 
splitting experiment cannot be used to tune the qubit 
and the resonator exactly in resonance. 

In general however, noise is not white. Though the 
ohmic model which leads to constant relaxation rates 
k(u) and 7(w) is usually valid, the transition rates 7I 
coming from the dephasing bath can be asymmetric. 
This yields a contribution r\$ to the total asymmetry 
V = ^us + f}<j>i where ?7 US is given by equation (391 and 



1-4A 2 



Vj> 



[7^(Ai_,i + )-7^(A 1+il _)] 



(40) 



As discussed in Appendix |B 2 1 noise at negative frequen- 
cies will only appear for non-zero effective bath tem- 
peratures, i.e. the rates respect detailed balance. With 
Ai-t^izp = ±2g, we therefore obtain 

70 (-2 3 ) - cxp (-j^) 70(25), (41) 



and have 



V<t> 



1-4A 2 



" M> ! -MT ' " ' 



We now distinguish two limits. If /c_bT ^> 2g, 77^ — > and 
we retrieve the classical noise limit. The asymmetry is 
then entirely due to the ultrastrong coupling. It is possi- 
ble to isolate this ultrastrong signature by increasing the 
effective temperature of the bath, for example in circuit 
QED, by injecting noise in u a with an external flux line. 
In the opposite scenario, if kgT <C 2g, the asymmetry 77^ 
becomes important. In particular, if T — > 



1 - 4A 2 .„ . 
V<t> - — 1A 2 9>- 



(43) 



70(2<?)- (42) 



Knowing the ultrastrong contribution to asymmetry, ei- 
ther by calculating it with equation (39) or by measuring 
it experimentally, it is possible to isolate the effect of 
quantum dephasing noise by taking 77^ = 77 — tj us . This 
asymmetry is thus a probe for dephasing noise at the vac- 
uum Rabi splitting frequency. The ultrastrong coupling 
regime widens the range of accessible values of g. As a 
result, the noise spectrum entering the rate 7</>(w) could 
realistically be investigated to frequencies up to ~ 2 GHz, 
where data is lacking [35] and where a crossover from 1// 
to ohmic behaviour is expected to happen |49] , 



B. Qubit frequency modulations: sidebands and 
photon generation 

In this section, we focus on the effect of the term f(t)a z 
of the classical dephasing model Eq. (13). We first con- 
sider the case where /(£) is a controlled modulation of 
the qubit frequency (for example, using an external flux) 
before turning to the situation where f(t) represents inco- 
herent noise. Both cases will be related to the dynamical 
Casimir effect [SU] [ST]. 

For g <C S, A, we apply on H = Hr + f(t)a z the 
dispersive transformation [52j . here generalized to take 
into account the counter-rotating terms [37] 

Uu = exp { A (aV- - aa + ) + A (aa- - cr + a f ) } . (44) 

To second order in g, we find 

H D (t) ~H' + X '(tWaa z + f(t)a z - 2f(t)(\I + + A/ CR ) 
-2/(t)AA ( r z (a 2 + a t2) ) (45) 

where x'(t) = -2(A 2 + A 2 )/(£) and H' = [u r + (x + 
/j,)o~ z ]a 1i a+ [cj Q +x + /i]cr 2 /2 the free but Lamb and Bloch- 
Siegert-shifted Hamiltonian. 

We first focus on the case of a classical modulation 
f(t) = e z cos ujdt of the qubit transition frequency. For 
e z <C u>d, the oscillating terms proportional to a z and 
a)aa z can be dropped under the RWA. Depending on 
the choice of modulation frequency uj^, it is possible to 
select different terms in Ho(t) while dropping others. 
First, for cud = A, we have Hrj(t) ~ H' — e z \I + cor- 
responding to a red sideband transition. For u>d = S, 
we rather find H D (t) ~ H' — Ae z JcR corresponding to 
a blue sideband transition. Interestingly, these sideband 
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FIG. 5. Photon generation due to dephasing with the Lind- 
bladian Eq. (16 1. Full line : white noise. Dotted line : white 



noise with a cut-off frequency increasing from bottom to top. 
For the bottom dotted line the cut-off is such that only tran- 
sitions up to |i,±) for i — 2 are driven. For the top curve, 
transitions from \g0) through i — 8 are driven. Inset : photon 
generation rate j3 as a function of g for white noise. Points : 
numerical results. Line : perturbation theory Eq. ( |46[ ). The 
parameters are uj a /2n = uj r /2ir — 6 GHz, g/2n — 1 GHz and 
70/2tt = 1 MHz. 



rates are in first rather than second order in g/u>d- This 
is to be contrasted to the usual second order results ob- 
tained in circuit QED 53-55] and could be used to speed 
up two-qubit gates. Finally, modulating at twice the res- 
onator frequency, u>d = 2cu r , the Hamiltonian reduces 
to Hjj(t) = Hq — e z XAa z (a^ 2 + a 2 ), corresponding to a 
pumped parametric oscillator |25j . Rather than modulat- 
ing the resonator frequency |56j , here it is the qubit that 
acts as a moving boundary condition. In Ref. [41] this 
Hamiltonian was also linked to the dynamical Casimir 
effect. 

We now move to the situation where /(£) is a random 
function representing a classical dephasing bath whose 
spectral content may contain one or more of the above- 
mentioned relevant frequencies. If the spectral content 
extends to very high frequencies, it may act on the sys- 
tem through a combination of the above blue and red 
sideband transitions and photon-pair production, bring- 
ing it to an excited state which may display some degree 
of squeezing. While this discussion is only valid in the 
dispersive regime, we can extend these results to arbi- 
trary ratios g/ A. Indeed, in general the a z bath can 
drive any transition between Rabi eigenstates of same 
parity, as illustrated in Fig. [2] If the corresponding fre- 
quencies are present in /(£), combinations of qubit and 
resonator excitations are produced. In the simplest case 
where relaxation is neglected and the dephasing bath is 
white, this leads to a photon creation rate /3. As shown 



in appendix [El to second order in g, 
/3 = 2 7 ^A 2 T(9 2 ), 



(46) 



where we have defined T(6 2 ) = 1 + 2 cos 2 9 2 sin 2 6 2 . This 
expression is compared to exact numerics in the inset of 
Fig. [5J It analytically explains the A 2 dependence of the 
photon creation rate observed numerically by Werlang et 
al. [23J for the special case of cu a = uj r . In this work, the 
authors have used the standard Lindladian Eq. ( fl5| ), as- 
suming white noise and corresponding to the full line in 
Fig.[5j If the noise causing dephasing has a 1// spectrum, 
the standard Lindbladian therefore greatly exaggerates 
this effect. If noise decreases at higher frequencies, pho- 
ton generation has a smaller rate, but should also sat- 
urate, as illustrated by the dotted lines in Fig. [5] This 
is again a clear demonstration of the breakdown of the 
standard approach to treating dissipation in the presence 
of the counter-rotating terms. 

Finally, since the master equation is exactly the same 
in the quantum treatment of dephasing shown in Ap- 



pendix B 2 the above results remain valid in that case. 
However, the quantum approach explicitly incorporates 
temperature in a way that respects detailed balance. 
This implies that, at T — 0, 70 (w) = for lu < 0. Since, 
as shown in Appendix [El this Casimir-like photon gen- 
eration needs negative frequencies, a quantum dephas- 
ing bath could not generate excitations in the system at 
T = 0. In this model, photon production through de- 
phasing is thus intrinsically a thermal effect. 



VI. CONCLUSION 

We have shown the importance of treating the qubit- 
resonator system as a whole when studying its interaction 
with the environment. In particular, we have shown that 
the description offered by the standard master equation 
can break down, for example producing spurious qubit 
flipping or photon generation, even at zero temperature. 
To cure these unphysical problems, we have included the 
qubit-resonator coupling in the derivation of the mas- 
ter equation. The rates entering the modified master 
equation then depend on the spectrum of noise evaluated 
at the dressed transition frequencies. These rates have 
been obtained analytically for a qubit-resonator coupling 
g that is large enough for individual qubit and/or res- 
onator transitions to be resolved, and for the dispersive 
(\uj a — uj r \ 3> g) to the Bloch-Siegcrt (uj a + u r ^> g) 
regime. Even when including the counter-rotating terms 
in the qubit-resonator coupling, the results obtained here 
can be used beyond these regimes by relying on simple 
numerical diagonalization of the Rabi Hamiltonian. Re- 
sults in the ultrastrong coupling regime {g ~ w a , uj r ) have 
been presented. 

In our model, noise that caused pure dephasing in the 
standard master equation can now cause transitions in 
the system. In this sense, the master equation devel- 
oped here can be viewed as an extension of the dressed- 



dephasing model [32J0H]. I n the Bloch-Siegert regime, 
we find that the vacuum Rabi splitting spectrum can be 
asymmetric. This asymmetry can be used as a probe 
of the dephasing noise spectral density at currently un- 
explored frequencies ~ 1 GHz and above. Additionally, 
modulations of the qubit transition frequency can be used 
to generate red and blue sidebands, or as a parametric 
oscillator inducing squeezing. Finally, while this means 
that noise in a z can generate photons [531 E3 EH] i our 
model reasonably shows that these spurious excitations 
cannot be generated at zero temperature. 
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Appendix A: Dissipators for X and a x baths 

In this appendix, we derive in the dressed basis the 
Lindbladian corresponding to coupling to the X and a x 
baths. We take an arbitrary qubit-resonator system with 
the only assumption that the total excitation number has 

we assume 



a well-defined parity. As stated in section III 



that the system is coupled to two independent baths of 
quantum harmonic oscillators with an interaction of the 
form given by Eq. (12 1. Focusing here on only one bath, 



we find that in the interaction picture with respect to the 
free system and bath Hamiltonians, the coupling takes 
the form 



■c^)e- iHst {bie- ivit + b\e iuit ). 



( Al ) 
Expressing the system Hamiltonian in the dressed basis 

H s =J2E j \j)Ul (A2) 

3 

we have 

H SB (t) =Y, ai C jk \j)(k\ (bte-^ + b\j vtt ) e lA ^ 



jkl 



(A3) 



where Cj k = (j\(c + c^)\k) and Aj k = Ej — E k . We now 
split the sum in three parts 



JTsbW =$>' c iilJ>01 (toe" 

1,3 



-ivit 



b\e iVlt 



+ 1EE+EE \^ 3k \j)(k\ 

l j,k>j I j,k<j ) 

x (bie- liui - A i k)t + b]e i ^ l+A ^ i ') . (A4) 



Since C k j = C* k , this becomes 
ff SB(*) =EE a < C ^WI (bie-^ + b, 



tgii^t 



3 I 



+ E E a ' c ^ij)< fc i^ e ^ +Ajfc)f+h - c - ( A5 ) 

3,k>j I 

We now introduce the operator n = (— l) a a + a + a - 1 
whose eigenvalues label the parity of the total excitation 
number in the qubit-resonator system. Since [Hr,TL] = 
0, system eigenstates \j) have a well-defined parity. As c 
and c> change the excitation number by one, c + c* flips 
the parity when applied on a state. Therefore, Cjj = 0, 
which simplifies the Hamiltonian to 



where 



H SB (t) = s(t)B\t) + «t(*)B(t), (A6) 

*(*)= E C 3*\J){k\e iA * kt (A7) 

B(t)^Y, a ^ e ~ lult - ( A8 ) 



This formulation makes it easy to write the Born master 
equation for the system. Indeed, following the standard 
procedure, we find [39 



p I (t)= / d« / [*(« , )Pz(* , )*(*)-a(t)*(* / )Pi(«')K5 t (t)Bt(f)) 

Jo 

+ f&t' [ s \t') Pl {t') s \t) s\t) s \t') Pl {t')\ (B(t)B(t')) 
Jo 

+ f At 1 [s\t') Pl {t')s{t) - 8(t)s\t')pi{t')] (B\t)B(t')) 





+ / dt' [ S (t') Pl (t')s\t) - S \t) S (t') Pl (t')} {B(t)BHt')) 
Jo 

+ h.c. (A9) 

From this point, we make assumptions that are standard 
in the Born-Markov treatment of dissipation [3^ , except 
for the following considerations. In each term of the Born 
master equation, we find oscillating exponentials of the 
form exp[i(Aj k — Aji k i)t]. Since k > j and k' > j', the 
argument of these exponentials will be zero for j = j' 
and k = k' ', or for pairs of different transitions in the 
system occuring at the same frequency. As discussed 
in Appendix [Cl in practice we are often interested only 
in a subset of the energy levels of the system for which 
all transitions have different frequencies. In that case, 
we can neglect all fast oscillating terms to obtain the 
following master equation in the Schrodinger picture 

p{t) = -i[H' s , P {t)] (A10) 

+ ]T P k n(A kj ,T)V[\k){j\])p(t) 

3,k>j 

+ J2 V k a + n(Ak 3 ,T))V[\j)(k\] P (t), 

3,k>j 
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with D*'* = 2ird(Akj)a 2 (Akj)\Cjk\ 2 and where we have 
introduced the density of states d(v) of the bath. We 
have also defined the Lamb-shifted system Hamiltonian 

H' s (t) = H s -J2\O j k\ 2 {L j k\k)(k\ 

+L' jk (\k)(k\-\j){j\)}- ( Al1 ) 

L': k are Lamb shifts caused by coupling to the environe- 
ment and are given by 



L 



jk 



L'; 



2tt J v + A jk 

P 



jk 



2tt 



dv 



T{v)n{v,T) 



(A12) 
(A13) 



v + Ajk 

and P is Cauchy's principal value. The function T(is) is 
a relaxation rate. In the case of photon loss, c —> a and 
we replace T(v) by n(v). In the case of qubit relaxation, 
c — > <r_ and we replace Y[v) by "f(v). In the main body 
of this paper, we neglect these Lamb shifts. 



Appendix B: Dissipators for the a z bath 
1. Classical model 

In this section, we derive the dephasing part of 
Eq. (16 1. For this, we introduce a stochastic function 



/(£) modulating the qubit frequency 

H dcp = f(t)<r z , 

where the mean value of /(£) vanishes. Following Ap- 
pendix [Al we express the Hamiltonian in the dressed ba- 
sis and move to the interaction picture with respect to 



Eq. (|Blj) to get 

H dcp (t)=f(t)Y / \j)(k\UW z \k)e^ t 

jk 



Following closely Ref. [3 
its Fourier decomposition 



/(*) 



(B2) 
we express f(t) in terms of 



dcu f(uj)e l 



to obtain 

H dcp (t) = Y / °l k \j)(k\f-A Jk (t), 
jk 

where we have defined 

/•Ajk + Bjk 
fA lk (t) = 



duj /(lj)^-^^. 



(B3) 

(B4) 

(B5) 



In writing this expression, we have considered that the 
main contribution to dephasing comes from a small fre- 
quency interval 2Bjf c around Ajk- For this approxima- 
tion to be valid, we must have Bjk <C Ajk- Using the 
Wiener-Khintchin theorem [59 



£[/(«)/(-«')] = *(« -«')£/(«), 



(B6) 



where E[x] is the classical mean value of x and Sf(ui) the 
spectral density of /(£), we then write 



/(w) = y^Mtfw), 



(B7) 



with £(v) such that E[£(w)] = and E[€(w)£(u/)] = 
5(u — u/), i.e. white noise. We now take Sf(ui) to be ap- 
proximately constant over each individual Bjk and con- 
sider that these bands do not overlap, allowing us to write 



f Ajh (t) = ^S f (A 



duj£(uj + A jk )e luJt . (B8) 



\jkj J ^ >L- --.I- I ->, 

-Bj k 



Assuming the dephasing timescale to be much slower 
than 1/Bjk, we can take Bjk — > oo and get 



U jk (t) = yJs f (A jk )U jk (t), 
finally leading to 



(B9) 



H dcp (t) =J2*i k \J){k\t-A dh (t)^S f (-A jk ). (BIO) 
jk 

If the transition frequencies Ajk are well-separated, we 
can treat each term of the above summation as an inde- 
pendent noise. This last form for H dop (t) then yields the 
following terms in the master equation 



( B1 ) J2h^' A ^ a ^ 2v ^^ + h^ v 



jk 



with 7^(-Aj-fc) = 2Sf(-A jk ). 



3 

(Bll) 



2. Quantum model 

To model pure dephasing in a quantum way, we intro- 
duce a quantum bath [39] 



H B =^vjb)bj 



(B12) 



The interaction of the system with this bath is given by 

H SB = ^2 <Xjkb]jbk<T+(T-, (B13) 

jk 

where bj is the ladder operator for bath mode j, with 
frequency Vj, and ajk is a coupling constant. This inter- 
action corresponds to the transfer of an energy quantum 
from one bath mode to the other through virtual excita- 
tion of the qubit. We now move to the dressed interaction 
picture with respect to Hg + H B 



Hj(t) =J2ajkb]b k e i(,/ i- /h)t e iHst <T-a + e 



iH s t 



(B14) 



jk 
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Using the closure relation of the system, we get 
H I (t)= Y, a 3^M i{Vi ~ Vk)t Zrnn\rn){n\e iA ™\ (B15) 

jkmn 

where we have defined the parity-preserving overlap ma- 
trix 



Z m „ = (m\a + <j-\n) 



(B16) 



In Appendix \K\ to obtain Eq. ( A10 ) for the coupling to 
the X and a x baths, we exploited the fact that all bath 
operators interacting with the system had zero mean 
value in thermal equilibrium |39j . This is not the case 
with the above interaction Hamiltonian, since terms for 
which j — k have a non-zero thermal mean value. To 
solve this problem, we include these terms in the sys- 
tem part of the total Hamiltonian, defining an effective 
shifted Hamiltonian 



where 

z mn (t) = Z mn \m){n\e iA ™ '. (B18) 

Assuming the bath is in thermal equilibrium, we get 

J3& = H s + J2 a J3 n {T)z mn {t). (B19) 



jmn 



We can now write the interaction Hamiltonian as 



with 



Hj(t) = B(t)s(t), (B20) 



B(t) = E a jk b]he i to- Vh)t (B21) 

S (t)=J2 Z mn(t). (B22) 

mn 

This allows to use the Born master equation 

pi(t) = - I dt' (B23) 

Jo 

x {[ s (t)s(t') Pl (t') - s(t') Pl (t') S (t)] (B(t)B(t')) p 
+ [ Pl (t') s (t')s(t) - s(t) Pl (t'W)] (B(t')B(t)) ) , 

where /3 indicates that the mean value is taken in a ther- 
mal state. The above correlators take the form 

(B(t)B(t')) = J2 a^Ji^il + MT))^-^ 

(B(t')B(t)) = J^ a^(T)(l + fT fe (T))e-^-^^, 

(B24) 



where r = t—t' and where we have taken the system-bath 
coupling constant to be real and symmetric under the 
exchange of modes j and k. The Born master equation 
becomes 

pi(t) = - E 

mnm'n' 
x [ [z mn (t)z m 'n'(t)pi(t - T) - Z m 'n'(t)Pl(t ~ r)z mn (t)) 

T (B(t)B(t-r)) 
+ [ P i(t-T)z m , n ,(t)z mn (t) - z mn (t)pi{t-r)z m 'n'(t)] 



x / dr e 
'o 



-iA„ 



: I dr e" 
/o 



-iA, 



'^(B(t-r)B(t))}, 



(B25) 



Replacing pj(t — r) by pi(t) and extending the upper 
boundary of the integrals over time to infinity, i.e. doing 
the Markov approximation, we get 



#S = #8 + l>ii&J&i*m»(*)> ( B17 ) / dTe " 

" 



dr e" 



with 

Tmn = 27T 

and 



*A m ,„,r ( B {t)B{t - T)) ~ ^7m'n'(T) - %L m , n , 
!A ™'"' 1 (B(t - T)B(t)) ~ \lm'n'(T) - iL m , n /, 

(B26) 

dv a 2 (v, u + A mn )d(u)d(v + A m „)x 
n(i/,T)(l + n(i/ + A mn ,T)), (B27) 



^ - *> - A m „ 

xn(v,T)(l + n(p',T)), (B28) 

As in Appendix [AJ we assume that all relevant tran- 
sitions in the system have different frequencies. This al- 
lows to drop fast rotating terms. Conditions in which 
the present approach might be inaccurate are explained 
in Appendix |C| 

Knowing that 



Z,, 



(B29) 



we obtain the following master equation in the 
Schrodinger picture 



p(t)=- l [Hl P (t)] + 1 -Av 



E<"»h 



Pit) 



E 



7* (A 



n.m) I m n\2 

\cr. 



V[\m){n\\p{t), (B30) 



with 70(A„ m ) = 7 nm /2, and where we have defined the 
Lamb-shifted Hamiltonian 



H's =H' S + E \ Z ™\ 2 Lran\n)(n\ 



(B31) 
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Equation ( B30 ) is exactly the master equation found for 



a classical bath if we neglect Lamb shifts, which can be 
done at low temperature and system-bath coupling. 

Finally, since the above master equation has been ob- 
tained for a bath in thermal equilibrium, the rates must 
obey detailed balance [47] 



7^(-w) = cxp - 



knT 



70 (w). 



(B32) 



Appendix C: Conditions under which the master 
equation developed here is applicable 

Here, we discuss conditions under which all relevant 
transitions have different frequencies and the above mas- 
ter equation can safely be applied. As stated in Sec- 
tion |IV A| in the dispersive regime, if £ ~ k, resonator 
transitions overlap. On the other hand, if the ratio g/A 
is large enough to have ( ^> k, this degeneracy is lifted, 
at least for low excitation numbers. Indeed, for high exci- 
tation numbers, some transitions might accidentally have 
the same frequency. 

We now define a critical excitation number n cr ; t under 
which all transitions occur at different frequencies and 
thus can be treated independently, given a sufficient ratio 
g/A. We limit ourselves to the Bloch-Siegert regime, 
under which g <C S and the counter-rotating terms are 
treated in a perturbative way. In this case, the energy 
levels are 



E n ± ~ nu> r ± 



1 



(ABS)2 + 4 5 2^ 



(CI) 



and thus display a nonlinearity scaling in yfn. In ad- 
dition, and as illustrated in Fig. [61 the bath operators 
only couple states in the same doublet, or one or two 
doublets away from each other. Moreover, as explained 
in Sec. |IIIB| parity selection rules apply so X and a x 
baths can drive only transitions between adjacent dou- 
blets (Type 1), while o~ z noise can induce transitions 
inside a doublet (Type 0) or between second-nearest- 
neighbour doublets (Type 2). This allows us to find a 
distinct n CT i t for individual baths by looking at every pos- 
sible combination of transitions and finding when some 
possibly overlap. 

We first consider transitions driven by X and o~ x baths, 
which are of Type 1. At low n, these essentially occur at 
distinct frequencies because of the yfn nonlinearity. How- 
ever, as n increases, y/n becomes comparable to y/n + 1 
and these transitions become closer. When the frequency 
difference between two transitions becomes of the order of 
nK, their typical linewidth, our model breaks down. We 
get an order of magnitude estimate of this critical n with 
the condition (E n+ i t+ — E n<+ ) — (E n +x- — E n ,—) ~ nn. 
Dropping terms of order higher than g 2 and assuming 
g ^$> k and g > |A|, i.e. such that the system is out of 




1000 
g (MHz) 



FIG. 6. a) Types of transitions allowed in the Bloch-Siegert 
regime (g <C E). Red: transitions driven by the (even) a z 
bath. Black: transitions induced by the (odd) X and o x 
baths, b) Exact energy levels of the Rabi Hamiltonian with 
increasing coupling strength, obtained numerically. Crossings 
between levels in the spectrum lead to pairs of transitions 
with equal frequency. However, since in each of these pairs, 
one transition is even and the other is odd, they belong to 
different baths and these overlaps are not relevant for the 
master equation. Parameters are io r /2-K = u> a /27r = 6 GHz. 



the dispersive regime, we get 



"crit 



A 
25 



-i 2/3 



(C2) 



Typically, this means we can have hundreds of excitations 
before the nonlinearity disappears and the model breaks 
down for odd baths. As a result, this limitation is not 
relevant in practice. 

We now turn to <7 z -driven transitions. Type 2 transi- 
tions will start to overlap in similar conditions as above, 
but with higher n, since they involve more widely sep- 
arated energy levels. Thus, these transitions do not 
set h cr it ■ Because of the y/n dependence of the non- 
linearity, all Type transitions have different frequen- 
cies. However, Type and 2 transitions can overlap. 
These overlaps depend on parameters g, uj r , u> a , and 
n in a nontrivial way, but a critical n for which this 
starts to be possible can be established with the condition 
E n +i,- — -En-1,4- = E n , + — E n ,-, which leads to a compli- 
cated expression for n cr it. Yet, we can get an estimate for 
this critical number with the criterion E n ^ + — E n _ ~ co r , 
which leads to 



'(0-2) 



4(g 2 + /lA) ' 



(C3) 



In resonance with uj a /2ir — u r /2ir — 6 GHz and for 
g/2-K = 1 GHz, rv cl ^ t — 9, enough to accurately de- 
scribe many experiments such as spectroscopy. 
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Finally, we emphasize that over n clit , we can only 
say that some pairs of transitions might overlap. If the 
involved levels do not play a leading role in the dynamics 
of the system under study, the master equation presented 
here should still give sensible results in practice. 



Appendix D: Transition matrix elements 

In this Appendix, we evaluate the overlap matrix ele- 
ment between eigenstates \j) and \k) of the Rabi Hamil- 
tonian for an arbitrary operator O. This is done using 
the perturbation theory presented in Section [TT] to second 
order in g such that 



O lk ~ {j\lfOU\k), 



(Dl) 



where \j) is the eigenstate in the Bloch-Siegert basis, 
Eq. ([9|. Using Eqs. ([6| and with the unitary trans- 
formation Eg . pi , we calculate the transition matrix el- 
ements Eqs. pM, d22| and d23l. For O = X, we obtain 



(D2) 



X 90;i- = (i + l)sme 1 -lcose 1 
X9°'A+ = (l + l)cosd 1 +lsmd 1 

X n+,n+l,+ = ryj^j _ q s[n9n + l cos n ] Sm 0„ +1 
+ Vn + 1(1 + I) COS n COS n +l 

X n+ > n+1 '- = - [y/n(l - I) sm9 n + I cos 9 n ] cos 9 n+1 



+ vn + 1(1 + I) cos 9 n sin n +i 
X n -- n+1 '+ = [-Vn(l - 0cos6»„ + 1 sin n ] smO n+1 

+ vn + 1(1 + ?) sin 9 n cos #„+i 
A" +;n+1 '- = - [- Vn(l - cos 6>„ + ? sin 9 n ] cos 0„ +1 

+ yn + 1(1 + I) sin$„ sin# n+ i, 

where I = 2£ + l 2 /2; £ = gA/2uj r is defined below Eq. ||| 
We note that Xjj = Xj, b . All other matrix elements are 
zero to second order. Similarly, for O — o~ x we find 



r 90;l- 



s sin &i 



(D3) 



af ;1 + = -r 2 , sin 6»i - s cos6»i 



ra— ;n+l,- 



[r 2 , sin 9 n+1 + s n+1 cos 9 n+1 ] cos 0„ 
[s„ sin 9 n+ i + t„ cos n+1 ] sin 9 n 

l] COS 6* 

s n cos 6< n+ i + £„ sin 9 n+ {\ sin 6> n 

Sn+iCOS^n + i] sill 6>„ 

- [s n sin n+1 + t n cos 0„ +1 ] cos 6>„ 

- [—r„ cos 6»„ + i + s n+ i sin n +i] sm # 

- [-s n cos9 n+ i +t n sm9 n+ i]cos9 n , 



_n+;n+l,+ 



< +;n+1, = - [-r£ cos (9 n+ i + s n+1 sin 



ar ; " +1 ' + = -Ksm^+i 



with r n = 1 - A 2 (n + 1/2), s„ = Ay 7 ™, and i„ = 
2^v /n ( n + !)• Finally, O = cr z yields 

af ;9 ° = 2A 2 - 1 (D4) 



CTf ;2+ = 2Asin6» 2 
crf ;2 " =-2Acos0 2 



n+;n+ 



[2A 2 (n - 1) - l] cos(2(9„) + 4A 2 cos 2 9 r , 



a n+;n- = 2 (2A 2 n - 1) sin 9 n cos9 n 



[2A 2 (n - 1) - l] cos(2(9„) + 4A 2 sin 2 9 n 



T n+;n+2,+ _ 



= 2Ayn + 1 cos 9 n sin 9 n +2 



cr" +;n+2 '- = -2AvVTT 



COS U„ COS 



n COS P n+2 



n-;n+2,+ _ 



= 2Ayn + 1 sin 9„ sin 0„+ 2 



n— ;n+2,— 



= — 2Ay/n + 1 sin 9 n cos 



n+2- 



Appendix E: Photon creation rate under white a z 
noise 

In this section, we derive the photon creation rate 



Eq. ( 46 ) caused by white noise fluctuations in the qubit 



transition frequency. To simplify the discussion, we con- 
sider only transitions to |2,±), the first accessible dou- 
blet. 

The photon creation rate is given by 



d 



fi = - (at a) = Tr [p(t)a*a] 



df 



(El) 



We take the initial state to be \g0). To obtain a constant 
rate, we limit ourselves to very small times i, such that 
j3 ~ /3(0). As illustrated in Fig. pH this will not cause 
any problem for white noise, for which (a)a) is found 
numerically to increase linearly at all times. Since p(0) = 
— i [Hr, p(0)] + £dr/?(0) and p(0) commutes with Hr, we 
obtain 



/3~Tr 



a^a£ dY \gO)(gO\ 



(E2) 



As shown in Eq. (16), £dr has a component responsible 



for pure dephasing and another that induces transitions. 
Since p(0) is an eigenstate, the dephasing term cancels 
out. We thus get, after moving to the Bloch-Siegert basis 



/3--(r0 2 -^ + r^ o )(. 9 o|( a t a )BS| ffO) 



(E3) 



r 



2-,s0/ 



|(at a ) BS |2-)+r^°(2+|( a t a )BS| 2+) . 



With 

{a^a) 35 = a? a - A(acr_ + a^cr + ) - 2£(a 2 + a t 2 ) CT 
1\ 1 



A z 



a)a+ - | n 



-A 2 , 
2 ' 



(E4) 



and using Eq. ( 19 1 for the transition rates as well as 



Eq. ( D4 1 for the corresponding transition matrix ele- 



ments, we obtain 

~ 2A 2 [T 2 _(0 2 ) 70 (- W2 _) + T 2+ (9 2 )^(-u 2+ )} , (E5) 
where 

T 2 _ = (1 + sin 2 9 2 ) cos 2 9 2 (E6) 
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T, 



2+ 



(l + cos 2 2 )sin 2 2 , 



(E7) 



and w 2 ± = E 2 ± — E, 



gO- 



Eq. ( E5 1 clearly shows that the spectrum at large neg- 



ative frequencies must be important in the a z bath for 
the photon generation rate to be significant. However, in 
our model, this bath respects detailed balance. Indeed, 
70 (— w) = exp(— uj/kBT)^^,(uj) such that 7^(— cj) — > 
for oj 3> ksT, meaning that these contributions should 
be very small for low temperatures. 

Finally, when 7(— u>2-) = t( — w 2+) = 7^j which is 
the case for white noise, Eq. |E5| reduces to the simpler 



Eq. (46) 



Appendix F: Vacuum Rabi splitting 



As outlined in Section V A here we calculate (a) . un- 



der weak cavity driving. We will assume that dephas- 
ing noise at high negative frequencies is weak, such that 
transitions from the ground state to the |2,±) doublet 
are negligible as shown in Appendix |EJ Together with 
the weak driving assumption, this means that only the 
first three levels of the system are relevant. 

For simplicity, we first move to the Bloch-Siegert basis 
defined by Eq. ^ to get 



#drvnW = U^H{t)U = H BS +ea BS e u " + h.c. 
It is also useful to move to a rotating frame with 

V(t) = e -^[( ata ) BS +^ s /2]* ) 
to obtain the time-independent Hamiltonian 



(Fl) 



(F2) 



Hl s vn = A B Va + ^fa z + gl + + e(a BS + a t BS ), (F3) 



with 



,BS 



i>r — ^d — \X 



A BS = io a - uj d + fi. (F4) 



The Heisenberg equation of motion for an arbitrary 
operator O is 



±(d)=i([H dTVD ,d 
where in the subspace {\g0), |1— ' 



C Ar O 



(F5) 



|1+)} 



£g- = ]T (n° M + Tf lCT ) V [|sO)(la 






ll-)(l+l 



-rl^-vo 



i+)(i- 



Vr 



+ $ 1 +> 1 +|i+)(i+| •, (F<n 



where V [Q]6 = (2 &6Q - &QO - 0&Q)/2. Rates 
are defined in Section lillBI 



We are interested in obtaining the mean value of a and 
cr_ . Since mean values does not depend on the frame, we 
will simplify calculations by working in the Bloch-Siegert 
picture. To do so, we first calculate the effect of the 
dissipators on a BS and ct5 s , knowing that 



BSiABS 



V [Q^]0 



* „\ BS 

Vo[Q]o) . 



(F7) 



This means that we can treat the states and operators 
in the effective Jaynes-Cummings Hamiltonian basis and 
then use the unitary U defined in Eq. pi) to move back 
to the Bloch-Siegert frame, which takes the non-RWA 
terms in consideration. In the three-level approximation, 
we have 

a= cos 0i 1^0) (1 + |+ sin 6>i |#0) (1-| (F8) 

o-_ =-sin(9 1 |. 9 O)(l + |+cos0 1 |. 9 O)(l-|, (F9) 

resulting in the dissipators 



Z>o[|flO><l- 

£o[|l-)(H 
^o[|l+)(l- 
V o [\g0){l + \]<r- 
T>o[\gO){l - \]<r- 
Z>o[|l-)<l + |]cr_ 
■D [\1+){1-\]t- 



- - cos 0i (a cos 0! — a_ sin 0i ) 

- - sin 0i (a sin 0! + er_ cos 0i ) 

- - cos 0i (a cos 0i — a- sin 0i ) 



sin 0i (a sin 0i + er_ cos 0i ) 
sin 0i (a cos 0i — cr_ sin 0i ) 
cos 0i (a sin 0i + c_ cos 0i ) 
sin 0i (a cos 0i — cr_ sin 0i ) 



cos 0i (a sin 0i + a_ cos 0i ) 



Proceeding in a similar way for dissipators involved in 
pure dephasing yields 



V% ph a = - 



\(4 



cos z 01 + ' 



sin f i J a 



2 
- - sin Vx cos f i l 7^ 

X>^ pl V_ = - - sin 0i cos 0i (7^ - jt) <i 



■ o 



2 V 7 * 
where we have defined 



2 ("ft sin2 "i 



+ 7^(0) \„ 9 0,g0 



cos f 1 a 



,i±,i±l 



(F10) 



(Fll) 



(F12) 



If we now add the contributions from all dissipators in 
Eq. (F5) to the Heisenberg equation for a, we get 



-r + (0i)(a)-r ? (0 1 )(cr_) 
while for cr_, we obtain 

-r ? (0 1 )(a)-r_(0i)(cr_) 



(F13) 
(F14) 
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Here, we have defined 



r+(0i) = Ti sin 2 9 Y + T 2 cos 2 B x 

T-(di) = r x cos 2 Ox + T 2 sin 2 Ox 

r)(8x) = (T 1 -r 2 ) sin 9 X cos 0x, 



with the rates 



t 



Ti 



7- + % + 7 



4- 



r? 



7+ + % + 



2 2 

This in turn involves the expressions 

7 ± =k(Ai ± , s o)|^°> 1± | 2 



(F15) 
(F16) 

(F17) 



(F18) 



^(Axi^kf^l 2 . 



vt/4. 



- 70 (A 1T , 1± )|^' 1± | 



(F19) 
(F20) 



We now calculate [-ffdrvn,a BS ] and [i/drvn, 0"? S ] ■ Ne- 
glecting terms that lead to leakage out of the effective 
Jaynes-Cummings three-level system, we obtain simple 
forms for a BS and ct bs in the bare basis 



,BS 



RS 



1 + — 1 a - Acr + + 2£a f 
1 + — J <t_ - Aa f . 
From these expressions, we easily get 

[Hdrvn,® ] — — e 



1 + ~W% 



1 



A 5 



(F21) 
(F22) 

3<7- 



BS\ 



(2£A^ - mK + (2& - AAf>+ (F23) 



[H dlvn ,a BS ] ~_ l 



A 5 



A BS a 



1 



A 2 



AA^a 1 - /io- 4 



,9a 



(F24) 



We now want to express this result in the Bloch-Siegert 
basis. From Eqs. (IF21J) and (|F22h 



1 



A 2 



RS 



AatBS. 



(F25) 



Knowing that a Bb =a + 0(A), Eq. (F21 1 leads to 



A 2 
1- T I» 



BS 



Act+ s - 2£a 



fBS 



(F26) 



This allows to find the appropriate commutators 



[-T^drvn; f* 



BS1 



- e - A BS a BS _ | j 



A 2 
2 



5^- 



2(£A BS 
[A(Af 



tBS 



A BS )- 2 e 5 K 



BS 



[F drvn , a BS ] ~ - [(1 + A 2 ) Af + A 2 A BS ] a 



[A(A B 



BS 



A 



BS\ 



2fr] 



+ BS 



,BSl „BS 

BS 
+ • 



2/1(7 



(F27) 

ga BS 

(F28) 



We can now write equations for the evolution of a 
and a BS . While Hamiltonian contributions are purely 
imaginary, those coming from dissipation are purely real. 
Imaginary terms lead to oscillatory behaviour, while real 
terms account for excitation and relaxation. Terms in 
a and a_ have both real and imaginary components, 
but terms in a^ BS and <r B only have imaginary contri- 
butions. The latter then only contribute through oscilla- 
tions in the dynamics. Since we are only interested in the 
steady-state behavior, we can neglect them. This allows 
to get the following equations for the steady state if we 
neglect terms of order higher than g 2 



iA BS + r+(0i)] (a BS ) + [ig + r,(6x)} (a BS ) +ie = 
[ig + V (9x)} (a BS ) + \iA BS + r_(0i)l (a BS ) = 



where we have defined A BS = (1 + A 2 )A BS + A 2 A BS . 
Solving the above set of equations, we finally obtain 



ieG 



siyij 



G»(0i) - G q (6x)G r (6xY 



(F29) 



where 



G q (9x) = r_(9 1 ) + iA BS 


(F30) 


G r (ex) = r + (9x) + iA BS 


(F31) 


Gr,(9x) = ig + v(0x). 


(F32) 



[1] S. Haroche and J.-M. Raimond, Exploring the Quantum: 

Atoms, Cavities, and Photons (Oxford University Press, 

Oxford, 2006). 
[2] M. Brune, E. Hagley, J. Dreyer, X. Maitre, A. Maali, 

C. Wunderlich, J. M. Raimond, and S. Haroche, Phys. 

Rev. Lett., 77, 4887 (1996)] 



[3] S. Gleyzes, S. Kuhr, C. Guerfin, J. Bernu, S. Delegfise, 

U. Hoff, M. Brune, J. Raimond, and S. Haroche, Nature, 

446, 297 (2007). 
[4] C. Guerlin, J. Bernu, S. Deleglise, C. Sayrin, S. Gleyzes, 

S. Kuhr, M. Brune, J. Raimond, and S. Haroche, Nature, 

448, 889 (2007). 



16 



[5] 
[6] 



181 

[9] 

[10] 

[11] 



[12] 
[13] 
[14] 

[16] 

[19] 
[20] 

[21 



A. Blais, R. S. Huang, A. Wallraff, S. M. Girvin, and 
R. J. Schoelkopf, Physical Review A, 69, 62320 (2004). 
A. Wallraff, D. Schuster, A. Blais, L. Frunzio, 
R. S. Huang, J. Majer, S. Kumar, S. Girvin, and 
R. Schoelkopf, Nature, 431, 162 (2004). 

A. A. Houck, J. A. Schreier, B. R. Johnson, J. M. Chow, 
J. Koch, J. M. Gambetta, D. I. Schuster, L. Frunzio, 
M. H. Devoret, S. M. Girvin, an d R. J. Schoelkopf, |PhysJ, 
Rev Lett., 101 08 0502 (2008)] 

~~T. Gambetta, A. Blais, D. I. Schuster, A. Wallraff, 
L. Frunzio, J. Majer, M. H. Devoret, S. M. Girvin, and 
R. J. Schoelkopf, Physical Review A, 74, 42318 (2006). 
J. Gambetta, A. Blais, M. Boissonneault, A. A. Houck, 
D. I. Schuster, and S. M. Girvin, Physical Review A, 77, 
12112 (2008). 

D. I. Schuster, A. Wallraff, A. Blais, L. Frunzio, R.-S. 
Huang, J. Majer, S. M. Girvin, and R. J. Schoelkopf, 
Phys. Rev. Lett., 94, 123602 (2005) 

F. R. Ong, M. Boissonneault, F. Mallet, A. Palacios- 
Laloy, A. Dewes, A. C. Doherty, A. Blais, P. Bertet, 

D. Vio n, and D. Esteve, |Phys. Rev. Lett., 106,T6 7002 
"(20TT)1 

J. Bourassa, J. M. Gambetta, A. A. Abdumalikov, 
O. Astafiev, Y. Nakamura, and A. Blais, Physical Re- 
view A, 80, 32109 (2009). 

T. Niemczyk, F. Deppe, H. Huebl, E. Menzel, F. Hocke, 
M. Schwarz, J. Garcia-Ripoll, D. Zueco, T. Hummer, 

E. Solano, et al, Nature Physics (2010). 

P. Forn-Diaz, J. Lisenfeld, D. Marcos, J. J. Garcia-Ripoll, 
E. Solano, C. J. P. M. Harmans , and J. E. Mooij,[Fhys] 
Rev. Lett., 105, 237001 (2010)| 

J. Casanova, G. Romero, I. Lizuain, J. J. Garcia-Ripoll, 
and E. Solano, |Phys. Rev. Lett., 105,1 63603 (2010)| 

B. Peropadre, P. Forn-Diaz, E. Solano, and J. J. Garcia- 



Ripoll, Phys. Rev. Lett., 105, 023601 (2010) 

J. Hausinger and M. Grifoni, Phys. Rev. A, 82, 062320 

(2010) 



J. Hau singer and M. Grifoni, Phys. Rev. A, 83, 030301 
T20lT)1 



J. Hausinger and M. Grifoni, New Journal of Physics, 10, 
115015 (2008). 



[22" 

[23; 

[24 

[25; 

[26] 
[27] 

[28] 
[29] 

[30] 

St 



C. Ciuti and I. Carusotto, Phys. Rev. A, 74, 033811 

(2006) 

S. De Liberato, D. Gerace, I. Carusotto, and C. Ciuti, 

|Phys. Rey. AT~8 0, 053810 (20097) 

P. Nataf and C. Ciutl^ ArXiv 

arXiv:1106.1159 



e-prints (2011) 



quant-ph] 



T. Werlang, A. V. Dodonov, E. I. Duzzioni, and C. J. 
Villas-Boas, |Phys. Rev. A, 78, 0538051 2008) 
M. O. Scully and M. S. Zubairy, Quantum Optics (Cam- 
bridge University Press, 1997). 
D. Walls and G. Milburn, Quantum Optics (Springer, 

Berlin, 1994). 

D. Braak, Phy s. Rev. Lett., 107, 100401 (2011)) 

I. Shavitt and L. Redmon, The Journal of Chemical 

Physics, 73, 5711 (1980). 



S. Ashhab and F. Nori, Phys. Rev. A, 81, 042311 (2010) 
M. Scala, B. Militello, A. Messina, S. Maniscalco, J. Piilo, 
and K. Suominen, Journal of Physics A: Mathematical 
and Theoretical, 40, 14527 (2007). 



R. Xu, Y. J. Yan, and X. Q. Li, Phys. Rev. A, 65, 023807 

(2002) 

M. Reed, B. Johnson, A. Houck, L. DiCarlo, J. Chow, 

D. Schuster, L. Frunzio, and R. Schoelkopf, Applied 



[32 



[33 



[34 



[35 

[36 

[37 

[38 



[39 
[40 

w 

[42 
[43 
[44 
[45 
[46 
[47 



[49; 
[50; 

[51 



Physics Letters, 96, 203110 (2010). 

A. Houck, D. Schuster, J. Gambetta, J. Schreier, B. John- 
son, J. Chow, L. Frunzio, J. Majer, M. Devoret, S. Girvin, 
et al., Nature, 449, 328 (2007). 



M. Boissonneault, J. M. Gambetta, and A. Blais, Phys. 
Rev. A, 77, 060305 (2008) 



M. Boissonneault, J. M. Gambetta, and A. Blais, Phys. 



Rev. A, 79, 013819 (2009) 

E. K. Irish, J. Gea-Banacloche, I. Martin, and K. C. 

Schwab, |Phys. Rev. B, 72, 195410 (2005)| 

E. K. Irish, [Phys. Rev. Lett 99, 1736 01 (2007)| 

D. Zueco, G. M. Reuther, S. Kohler, and P. Hanggi, 

Phys. Rev. A, 80, 033846 (2009) 



M.-J. Hwang and M.-S. Choi, Phys. Rev. A, 82, 025802 



(2010) 



H. Carmichael, An Open Systems Approach to Quantum 
Optics (Springer- Verlag, 1993). 



M. Boissonneault, J. M. Gambetta, and A. Blais, Phys. 
Rev. Lett., 10 5, 100504 (2010)] 

A. Dodonov, in Journal of Physics: Conference Series, 
Vol. 161 (IOP Publishing, 2009) p. 012029. 
M. Scala, B. Militello, A. Messina, J. Piilo, and S. Man- 
iscalco, |Phys. Rev. A, 75, 013811 (2007) | 

E. M. Purcell, H. C . Torrey, and R. V. Pound, |Phys. 
Rev., 69, 37 (1946) | 

F. Yoshihara, K. Harrabi, A. O. Niskanen, Y. Nakamura, 
and J. S. Tsai, |Phys. Rev. Lett., 97TT6 7001 (2006)) 

J. Fink, M. Goppl, M. Baur, R. Bianchetti, P. Leek, 

A. Blais, and A. Wallraff, Nature, 454, 315 (2008). 

X. Cao, J. Q. You, H. Zheng, and F. Nori, ArXiv e-prints 

(2010), arXiv:1009.4366 [quant-ph] 

A. Clerk, M. Devoret, S. Girvin, F. Marquardt, and 

R. Schoelkopf, Reviews of Modern Physics, 82, 1155 

(2010). 

J. Bylander, S. Gustavsson, F. Yan, F. Yoshihara, 

K. Harrabi, G. Fitch, D. Cory, Y. Nakamura, J. Tsai, 

and W. Oliver, Nature Physics (2011). 

A. Shnirman, G. Schon, I. Martin, and Y. Makhlin, 

|Phys. Rev. Lett., 9 4, 127002 (2005) 

J. R. Johansson, G. Johansson, C. 



M. 



and 



Wilson, 

F. Nori, Phys. Rev. Lett., 10 3, 147003 (2009) 
C. M. Wilson, T. Duty, M. SandbergTF. Persson, 
V. Shumeiko, and P. Delsing, Phys. Rev. Lett., 105, 



233907 (2010) 
[52] S. Haroche, J. Dalibard, J. Raimond, and J. Zinn- Justin, 

Fundamental Systems in Quantum Optics (Elsevier, New 

York, 1992). 
[53] A. Blais, J. Gambetta, A. Wallraff, D. I. Schuster, S. M. 

Girvin, M. H. Devoret, and R. J. Schoelkopf, Physical 

Review A, 75, 32329 (2007). 
[54] A. Wallraff, D. I. Schuster, A. Blais, J. M. Gambetta, 

J. Schreier, L. Frunzio, M. H. Devoret, S. M. Girvin, and 

R. J. Schoelkopf, Phys. Rev. Lett., 99, 50501 (2007). 
[55] P. J. Leek, S. Filipp, P. Maurer, M. Baur, R. Bianchetti, 

J. M. Fink, M. Goppl, L. Steffen, and A. Wallraff, Phys. 

Rev. B, 79, 1 80511 (2009) 
[56] T. Yamamoto, K. Inomata, M. Watanabe, K. Matsuba, 

T. Miyazaki, W. D. Oliver, Y. Nakamura, and J. S. Tsai, 

|Applied Physics Letters, 93 , 042510 (2008) 
[57] A. Dodonov, Physica Scripta, 82, 038102 (2010). 
[58] A. V. Dodonov, Journal of Physics: Conference Series, 

274, 012137 (2011). 
[59] C. W. Gardiner, Handbook of Stochastic Methods : for 

Physics, Chemestry and the Natural Sciences (Springer, 



17 
Berlin, 2004). 



